﻿using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Linq;
using System.Text;
using Common;

namespace ThreeItems
{
    public class DecentralizedAlg
    {
        /// <summary>
        /// 
        /// </summary>
        public double Tolerences { get; set; }

        /// <summary>
        /// 
        /// </summary>
        public double RepeatTimes { get; set; }

        /// <summary>
        /// 
        /// </summary>
        public double MaxIteration { get; set; }

        public void RunComputing(ExampleInfo example)
        {
            example.MaxQD1 = NormDistHelper.QNorm(example.UO1, example.E1, example.SD1);
            example.MaxQD2 = NormDistHelper.QNorm(example.UO2, example.E2, example.SD2);
            example.MaxQD3 = NormDistHelper.QNorm(example.UO3, example.E3, example.SD3);


            example.QD1 = NormDistHelper.QNorm(example.U1 / (example.U1 + example.O1), example.E1, example.SD1);
            example.QD2 = NormDistHelper.QNorm(example.U2 / (example.U2 + example.O2), example.E2, example.SD2);
            example.QD3 = NormDistHelper.QNorm(example.U3 / (example.U3 + example.O3), example.E3, example.SD3);

            //选取初始节点Q = 0
            //迭代
            bool stop = false;
            //存储连续小于容忍值的迭代
            int repeates = 0;

            for (int j = 0; ; j++)
            {
                if (stop)
                {
                    Trace.WriteLine("迭代结束!");
                    break;
                }

                if (j >= MaxIteration)//超出最大迭代次数
                {
                    Trace.WriteLine("超出最大迭代次数,结束!");
                    break;
                }

                double error = IterateOnce(example);
                if (error < Tolerences)
                    repeates += 1;
                else//出现任意一次不能接受，重置计数器
                    repeates = 0;

                //判断是否连续N次可以接受
                if (repeates >= RepeatTimes)
                    stop = true;

                //Trace.WriteLine(string.Format("第{0}次迭代,Error={1},QD1={2},QD2={3},QD3={4}", j, error, qd1, qd2, qd3));
            }

            Profit(example);

            Trace.WriteLine(example.PD);
        }

        public void Profit(ExampleInfo example)
        {
            //已经推导出来利润

            double variance1 = Math.Pow(example.SD1, 2);
            double variance2 = Math.Pow(example.SD2, 2);
            double variance3 = Math.Pow(example.SD3, 2);

            //计算有效需求的期望和标准差
            double pnorm1 = NormDistHelper.PNorm(example.QD1, example.E1, example.SD1);
            double gama1 = pnorm1 + (example.QD1 * (1 - pnorm1) - (variance1) * NormDistHelper.DNorm(example.QD1, example.E1, example.SD1)) / example.E1;

            double pnorm2 = NormDistHelper.PNorm(example.QD2, example.E2, example.SD2);
            double gama2 = pnorm2 + (example.QD2 * (1 - pnorm2) - variance2 * NormDistHelper.DNorm(example.QD2, example.E2, example.SD2)) / example.E2;

            double pnorm3 = NormDistHelper.PNorm(example.QD3, example.E3, example.SD3);
            double gama3 = pnorm3 + (example.QD3 * (1 - pnorm3) - variance3 * NormDistHelper.DNorm(example.QD3, example.E3, example.SD3)) / example.E3;


            double e1_e = example.E1 - example.R21 * example.E2 * (1 - gama2) - example.R31 * example.E3 * (1 - gama3);//有效需要的期望
            double sd1_e = Math.Sqrt(variance1 + Math.Pow(example.R21 * example.SD2 * (1 - gama2), 2) + Math.Pow(example.R31 * example.SD3 * (1 - gama3), 2));//有效需要的标准差

            double e2_e = example.E2 - example.R12 * example.E1 * (1 - gama2) - example.R32 * example.E3 * (1 - gama3);//有效需要的期望
            double sd2_e = Math.Sqrt(variance2 + Math.Pow(example.R12 * example.SD1 * (1 - gama1), 2) + Math.Pow(example.R32 * example.SD3 * (1 - gama3), 2));//有效需要的标准差

            double e3_e = example.E3 - example.R13 * example.E1 * (1 - gama1) - example.R23 * example.E2 * (1 - gama2);//有效需要的期望
            double sd3_e = Math.Sqrt(variance3 + Math.Pow(example.R13 * example.SD1 * (1 - gama1), 2) + Math.Pow(example.R23 * example.SD2 * (1 - gama2), 2));//有效需要的标准差


            double u1 = example.U1;
            double o1 = example.O1;
            double e1 = e1_e;
            double sd1 = sd1_e;
            double q1 = example.QD1;


            example.PD1 = (u1 + o1) *
                          ((e1 - q1) * NormDistHelper.PNorm(q1, e1, sd1) - sd1 * sd1 * NormDistHelper.DNorm(q1, e1, sd1)) +
                          u1 * q1;

            double u2 = example.U2;
            double o2 = example.O2;
            double e2 = e2_e;
            double sd2 = sd2_e;
            double q2 = example.QD2;
            example.PD2 = (u2 + o2) *
                          ((e2 - q2) * NormDistHelper.PNorm(q2, e2, sd2) - sd2 * sd2 * NormDistHelper.DNorm(q2, e2, sd2)) +
                          u2 * q2;

            double u3 = example.U3;
            double o3 = example.O3;
            double e3 = e3_e;
            double sd3 = sd3_e;
            double q3 = example.QD3;
            example.PD3 = (u3 + o3) *
                          ((e3 - q3) * NormDistHelper.PNorm(q3, e3, sd3) - sd3 * sd3 * NormDistHelper.DNorm(q3, e3, sd3)) +
                          u3 * q3;


            example.PD = example.PD1 + example.PD2 + example.PD3;
        }

        /// <summary>
        /// 一次迭代
        /// </summary>
        /// <param name="example"></param>
        /// <returns>误差e^2</returns>
        private double IterateOnce(ExampleInfo example)
        {
            double variance1 = Math.Pow(example.SD1, 2);
            double variance2 = Math.Pow(example.SD2, 2);
            double variance3 = Math.Pow(example.SD3, 2);

            //计算有效需求的期望和标准差
            double pnorm1 = NormDistHelper.PNorm(example.QD1, example.E1, example.SD1);
            double gama1 = pnorm1 + (example.QD1 * (1 - pnorm1) - (variance1) * NormDistHelper.DNorm(example.QD1, example.E1, example.SD1)) / example.E1;

            double pnorm2 = NormDistHelper.PNorm(example.QD2, example.E2, example.SD2);
            double gama2 = pnorm2 + (example.QD2 * (1 - pnorm2) - variance2 * NormDistHelper.DNorm(example.QD2, example.E2, example.SD2)) / example.E2;

            double pnorm3 = NormDistHelper.PNorm(example.QD3, example.E3, example.SD3);
            double gama3 = pnorm3 + (example.QD3 * (1 - pnorm3) - variance3 * NormDistHelper.DNorm(example.QD3, example.E3, example.SD3)) / example.E3;


            double e1_e = example.E1 - example.R21 * example.E2 * (1 - gama2) - example.R31 * example.E3 * (1 - gama3);//有效需要的期望
            double sd1_e = Math.Sqrt(variance1 + Math.Pow(example.R21 * example.SD2 * (1 - gama2), 2) + Math.Pow(example.R31 * example.SD3 * (1 - gama3), 2));//有效需要的标准差

            double e2_e = example.E2 - example.R12 * example.E1 * (1 - gama2) - example.R32 * example.E3 * (1 - gama3);//有效需要的期望
            double sd2_e = Math.Sqrt(variance2 + Math.Pow(example.R12 * example.SD1 * (1 - gama1), 2) + Math.Pow(example.R32 * example.SD3 * (1 - gama3), 2));//有效需要的标准差

            double e3_e = example.E3 - example.R13 * example.E1 * (1 - gama1) - example.R23 * example.E2 * (1 - gama2);//有效需要的期望
            double sd3_e = Math.Sqrt(variance3 + Math.Pow(example.R13 * example.SD1 * (1 - gama1), 2) + Math.Pow(example.R23 * example.SD2 * (1 - gama2), 2));//有效需要的标准差


            double nextQd1 = NormDistHelper.QNorm(example.U1 / (example.U1 + example.O1), e1_e, sd1_e);
            double nextQd2 = NormDistHelper.QNorm(example.U2 / (example.U2 + example.O2), e2_e, sd2_e);
            double nextQd3 = NormDistHelper.QNorm(example.U3 / (example.U3 + example.O3), e3_e, sd3_e);

            //判断误差Error

            //计算出来的概率
            double prob1_e = NormDistHelper.PNorm(example.QD1, e1_e, sd1_e);
            double prob2_e = NormDistHelper.PNorm(example.QD2, e2_e, sd2_e);
            double prob3_e = NormDistHelper.PNorm(example.QD3, e3_e, sd3_e);
            //待对比的概率
            double prob1_target = example.U1 / (example.U1 + example.O1);
            double prob2_target = example.U2 / (example.U2 + example.O2);
            double prob3_target = example.U3 / (example.U3 + example.O3);

            example.Prob1 = prob1_e;
            example.Prob2 = prob2_e;
            example.Prob3 = prob3_e;

            example.QD1 = nextQd1;
            example.QD2 = nextQd2;
            example.QD3 = nextQd3;


            ////已经推导出来利润

            //double u1 = example.U1;
            //double o1 = example.O1;
            //double e1 = e1_e;
            //double sd1 = sd1_e;
            //double q1 = example.QD1;


            //example.PD1 = (u1 + o1) *
            //              ((e1 - q1) * NormDistHelper.PNorm(q1, e1, sd1) - sd1 * sd1 * NormDistHelper.DNorm(q1, e1, sd1)) +
            //              u1 * q1;

            //double u2 = example.U2;
            //double o2 = example.O2;
            //double e2 = e2_e;
            //double sd2 = sd2_e;
            //double q2 = example.QD2;
            //example.PD2 = (u2 + o2) *
            //              ((e2 - q2) * NormDistHelper.PNorm(q2, e2, sd2) - sd2 * sd2 * NormDistHelper.DNorm(q2, e2, sd2)) +
            //              u2 * q2;

            //double u3 = example.U3;
            //double o3 = example.O3;
            //double e3 = e3_e;
            //double sd3 = sd3_e;
            //double q3 = example.QD3;
            //example.PD3 = (u3 + o3) *
            //              ((e3 - q3) * NormDistHelper.PNorm(q3, e3, sd3) - sd3 * sd3 * NormDistHelper.DNorm(q3, e3, sd3)) +
            //              u3 * q3;


            //example.PD = example.PD1 + example.PD2 + example.PD3;


            return Math.Pow(Math.Abs(prob1_e - prob1_target) / prob1_target, 2) +
                   Math.Pow(Math.Abs(prob2_e - prob2_target) / prob2_target, 2) +
                   Math.Pow(Math.Abs(prob3_e - prob3_target) / prob3_target, 2);
        }
    }
}
